# coding:utf-8
import itertools
import math
import numpy as np
import json
import csv

centers = [
    (139382.3437, 146047.0951),
    (157688.4185, 148489.0310),
    (149441.1622, 167278.6394),
    (131764.7130, 162103.0925)
]

traffic_block = [4.48, 6.57, 6.79, 6.53, 5.33, 6.35, 6.64, 6.78, 6.41, 6.79, 7.03, 6.96, 6.53, 6.08, 6.42, 6.64, 7.22,
                 6.49, 6.5, 5.4, 7.23, 6.33, 6.71, 8.01, 8.13, 6.28, 7.44, 7.71, 8.18, 9.21, 7.27, 6.29, 7.21, 6.78,
                 6.26, 6.32, 6.87, 6.09, 6.1, 6.95, 7.63, 6.99, 6.9, 6.72, 6.71, 6.69, 6.61, 6.13, 7.76, 6.12, 7.28,
                 6.72, 6.62, 7.66, 9.14, 10.94, 7.14, 9.21, 7.36, 8.03, 10.2, 6.6, 7.07, 6.03, 6.45, 6.41, 6.39, 6.32,
                 6.32, 6.41, 7.28, 7.36, 6.68, 6.62, 5.93, 2.93, 6.93, 6.27, 6.32, 6.42, 5.28, 6.26, 6.23, 4.06, 7.48,
                 7.39, 8.01, 8.14, 5.92, 8.69, 6.46, 7.07, 6.24, 6.26, 6.68, 6.73, 7.27, 11.54, 8.7, 6.04, 6.08, 6.17,
                 6.15, 5.94, 1.33, 3.91, 4.92, 5.94, 6.04, 5.94]

all_points = []

with open("points.txt", 'r') as f:
    lines = f.readlines()
    for line in lines:
        x, y = line.split(',')
        x, y = float(x), float(y)
        all_points.append((x, y))

data = json.load(open('out.temp'))
vs = data['fln']  # 点集
vs.sort()
out_mat = data['mat']  # 一级节点收发流量表
vs_pts = data['fln2pt']  # 一级节点对应的区域中心
vs_pts = [[pt - 4 for pt in pts] for pts in vs_pts]
pt2pt = data['edges']  # 区域中心收发流量表


def increase_by_year(year, rate=1.05):
    for i in range(len(out_mat)):
        for j in range(len(out_mat[i])):
            out_mat[i][j] *= (rate ** year)

    for i in range(len(pt2pt)):
        for j in range(len(pt2pt[i])):
            pt2pt[i][j] *= (rate ** year)

    for b in traffic_block:
        b *= (rate ** year)


# increase_by_year(30)

pt_ins = np.sum(pt2pt, axis=0)[4:]  # 区域中心流入流量
pt_outs = np.sum(pt2pt, axis=1)[4:]  # 区域中心发出流量

pt_ins_up = pt_ins * (4 / np.array(traffic_block))  # 区域中心可在地面上的最大流入流量
pt_outs_up = pt_outs * (4 / np.array(traffic_block))  # 区域中心可在地面上的最大发出流量

self_ins_up = [sum([pt_ins_up[pt] for pt in pts]) for pts in vs_pts]  # 一级节点可在地面上的最大流入流量
self_outs_up = [sum([pt_outs_up[pt] for pt in pts]) for pts in vs_pts]  # 一级节点可在地面上的最大发出流量

self_ins = np.sum(out_mat, axis=0)  # 一级节点的流入流量
self_outs = np.sum(out_mat, axis=1)  # 已经节点的发出流量


def get_distance(v1, v2):
    # 计算两个点的欧氏距离
    return math.sqrt((v1[0] - v2[0]) ** 2 + (v1[1] - v2[1]) ** 2)


def get_pairs(indics):
    # 计算点对组合（无重复）
    # param:indics:点的下标列表
    # return:pairs:点对
    pairs = list(itertools.product(indics, indics))
    pairs = filter(lambda x: x[0] < x[1], pairs)
    return pairs


def get_edges(vs, pairs):
    # 计算边
    # param:vs:点集;pairs:点对
    # return:edges[(a, b, d)]边集 a,b:点下标， d:距离
    edges = []
    for (a, b) in pairs:
        v1, v2 = vs[a], vs[b]
        distance = get_distance(v1, v2)
        edges.append((a, b, distance))
    return edges


def get_mst(vs, edges):
    # 生成最小生成树
    # param: vs:点集，edges:边集
    # return: mst:最小生成树边集
    ps = range(len(vs))
    mst = []
    for a, b, d in edges:
        p1, p2 = ps[a], ps[b]
        if p1 != p2:
            mst.append((a, b, d))
            for i in range(len(ps)):
                if ps[i] == p2:
                    ps[i] = p1
    return mst


def add_center_edges(mst, vs, centers):
    # 在最小生成树的基础上加入与物流中心之间的边
    # param: mst:最小生成树 vs:点集 centers:物流中心的点集
    # return: mst:加入与物流中心之间的边之后的边集合（不再是树）
    for i in range(len(vs)):
        tmp = []
        for j in range(len(centers)):
            tmp.append((i, len(vs) + j, get_distance(vs[i], centers[j])))
        min_edge = min(tmp, key=lambda x: x[-1])
        mst.append(min_edge)
    return mst


def get_links(mst, num):
    # 生成邻接表
    # param: num:节点数
    # return: links:邻接表
    links = [set([]) for i in range(num)]
    for a, b, _ in mst:
        links[a].add(b)
        links[b].add(a)
    return links


def get_paths(links, in_num, a, b):
    # DFS获取两点之间所有路径(有向)
    def get_paths_helper(links, in_num, a, b, path, paths):
        path.append(a)
        if len(path) > 1 and a >= in_num and a != b:
            return
        if a == b:
            paths.add(tuple(path))
            return
        for x in links[a]:
            if x not in path:
                get_paths_helper(links, in_num, x, b, path[:], paths)

        return paths

    return list(get_paths_helper(links, in_num, a, b, [], set([])))


def get_volumn(links, out_mat, num, in_num):
    # 获取每条边的流量
    # param: links:邻接表 out_mat:各个节点发出货物的关系 num:总节点数 in_num:一级节点数
    # return: volumn: 字典 key为边 value为流量
    volumn = {}

    for a, b in itertools.product(range(num), range(num)):
        if a != b and out_mat[a][b] > 0:
            paths = get_paths(links, in_num, a, b)
            path = min(paths, key=lambda x: len(x))
            for i in range(len(path) - 1):
                pair = (path[i], path[i + 1])
                volumn.setdefault(pair, 0)
                volumn[pair] += out_mat[a][b]
    return volumn


def get_cap(key, in_num):
    # 获取一条边的流量上限
    # param: key:边 in_num:一级节点数
    # return: cap:流量上限
    cap = 7200
    if any([n >= in_num for n in key]):
        cap = 14400
    return cap


def get_min_cap(key, in_num):
    # 获取一条边的最小流量上限
    # param: key:边 in_num:一级节点数
    # return: cap:流量上限
    cap = 3600
    if any([n >= in_num for n in key]):
        cap = 7200
    return cap


def get_edge_net(edge_volumn, in_num):
    # 建立边网络
    # param: edge_volumn:边流量字典 in_num:一级节点数
    # return: edge_net: 字典 key:边 value:(流量,流量上限)
    edge_net = {}
    for key in edge_volumn.keys():
        cap = get_cap(key, in_num)
        edge_net[key] = [edge_volumn[key], cap]
    return edge_net


def get_res_net(path, edge_net):
    # 建立残值网络
    # param: path:路径 edge_net:边网
    # return: res_net: 字典 key:边 value:(流量上限-流量,流量)
    res_net = {}
    for i in range(len(path) - 1):
        st = (path[i], path[i + 1])
        vol, cap = edge_net[st]
        res_net[st] = (cap - vol, vol)
    return res_net


def shunt_overflow_edge(a, b, links, edge_net, all_ins, all_outs, in_num):
    paths = get_paths(links, in_num, a, b)
    overflow = edge_net[(a, b)][0] - edge_net[(a, b)][1]
    if overflow < 0:
        return True
    for path in paths:
        forward = [edge_net[key][1] - edge_net[key][0] for key in zip(path, path[1:])]
        vol = min(forward)
        if edge_net[(a, b)][0] > vol:
            edge_net[(a, b)][0] -= vol
            for k in range(len(path) - 1):
                edge_net[(path[k], path[k + 1])][0] += vol
            for v in path[1: -1]:
                all_ins[v] += vol
                all_outs[v] += vol
            overflow -= vol
            if overflow <= 0:
                return True
    return False


def get_out_volumn(edge_net, i):
    # 计算某个节点的发出流量
    ans = 0
    for (a, b) in edge_net.keys():
        if a == i:
            ans += edge_net[(a, b)][0]
    return ans


def get_out_volumns(edge_net, num):
    # 计算所有节点的发出流量表
    ans = []
    for i in range(num):
        ans.append(get_out_volumn(edge_net, i))
    return ans


def get_in_volumn(edge_net, i):
    # 计算某个节点的流入流量
    ans = 0
    for (a, b) in edge_net.keys():
        if b == i:
            ans += edge_net[(a, b)][0]
    return ans


def get_in_volumns(edge_net, num):
    # 计算所有节点的流入流量表
    ans = []
    for i in range(num):
        ans.append(get_in_volumn(edge_net, i))
    return ans


def check(all_ins, all_outs, edge_net, in_num):
    vc_num = [0] * in_num
    for a, b in edge_net.keys():
        if edge_net[(a, b)][0] < 0:
            print "%d, %d negative vol" % (a, b)
            return False
    for a, b in edge_net.keys():
        if a < in_num and b >= in_num:
            vc_num[a] += 1
    for i, n in enumerate(vc_num):
        if n == 0:
            print "v %d has no center edge" % i
            return False
    for a, b in edge_net.keys():
        if edge_net[(a, b)][0] > edge_net[(a, b)][1]:
            print "(%d, %d) overflow, vol: %f, cap:%f" % (a, b, edge_net[(a, b)][0], edge_net[(a, b)][1])
            return False
    # 检测边网是否流量平衡
    for i in range(len(all_ins)):
        all_in, all_out = all_ins[i], all_outs[i]
        self_in, self_out = self_ins[i], self_outs[i]
        if abs((all_out - self_out) - (all_in - self_in)) > 1e-3 and \
                        all_out < self_out and all_in < self_in:
            print("vertex %d false" % i)
            print all_out, all_in, self_in, self_out
            return False
    return True


def calc_transport_rate(edge_net, in_num):
    # 计算转运率
    all_centers_in = [0] * in_num
    all_self_in = [0] * in_num
    for (a, b) in edge_net.keys():
        if a >= in_num and b < in_num:
            all_centers_in[b] += edge_net[(a, b)][0]
            all_self_in[b] += out_mat[a][b]

    transport_rate = []
    for i in range(in_num):
        rate = 0
        if all_centers_in[i] > 0:
            rate = 1.0 - all_self_in[i] / all_centers_in[i]
        transport_rate.append(max([rate, 0]))

    return transport_rate


def construct_net():
    # 初步构建物流网络

    # 生成最小生成树
    pairs = get_pairs(range(len(vs)))
    edges = get_edges(vs, pairs)
    edges.sort(key=lambda x: x[-1])
    mst = get_mst(vs, edges)

    # 添加与物流中心连接的边
    all_edges = add_center_edges(mst[:], vs, centers)

    all_edges.sort()
    # 计算邻接表
    links = get_links(all_edges, len(vs) + len(centers))
    # 计算边流量
    edge_volumn = get_volumn(links, out_mat, len(vs) + len(centers), len(vs))
    # 计算边网
    edge_net = get_edge_net(edge_volumn, len(vs))

    # 添加一级节点与物流中心的边
    for i in range(len(vs)):
        for j in range(len(centers)):
            edges.append((i, len(vs) + j, get_distance(vs[i], centers[j])))

    edges.sort(key=lambda x: x[-1])

    all_ins = get_in_volumns(edge_net, len(vs))
    all_outs = get_out_volumns(edge_net, len(vs))

    c = 0
    # 分流，直至所有的边都不会超出上限
    while any([edge_net[k][0] > edge_net[k][1] for k in edge_net]):
        print c
        keys = [key for key in edge_net.keys() if edge_net[key][0] > edge_net[key][1]]
        # 对于超出流量上限的边，不断根据残值网络进行分流，即将超出的流量分流到其他路径上
        for key in keys:
            print key
            # 按照从短到长的顺序将新边加入网络
            for a, b, d in edges:
                links[a].add(b)
                links[b].add(a)
                edge_net.setdefault((a, b), [0, get_cap((a, b), len(vs))])
                edge_net.setdefault((b, a), [0, get_cap((b, a), len(vs))])
                flag = shunt_overflow_edge(key[0], key[1], links, edge_net, all_ins, all_outs, len(vs))
                if edge_net[(a, b)][0] == 0 and edge_net[(b, a)][0] == 0:
                    # 如果这次分流中没有分流到新的边上，删除新的边
                    edge_net.pop((a, b))
                    edge_net.pop((b, a))
                    links[a].remove(b)
                    links[b].remove(a)
                if flag:
                    # 如果不在超出上限，停止迭代
                    break
        c += 1

    # 删除流量为0的单向边
    for (a, b) in edge_net.keys():
        if edge_net[(a, b)][0] <= 0.0:
            edge_net.pop((a, b))
            links[a].remove(b)
    # 添加单向边，保证都是双向边
    for a, b in edge_net.keys():
        if (b, a) not in edge_net:
            edge_net.setdefault((b, a), [0, get_cap((b, a), len(vs))])
            links[b].add(a)

    # 检测网络正确性
    print check(all_ins, all_outs, edge_net, len(vs))

    # 计算转运率
    transport_rate = calc_transport_rate(edge_net, len(vs))

    return links, edge_net, all_ins, all_outs, transport_rate


def get_neg_net(edge_net):
    neg_net = {}
    for a, b in edge_net.keys():
        neg_net[(a, b)] = edge_net[(a, b)][0]
    return neg_net


def reduce_edge_cap(links, edge_net, all_ins, all_outs):
    # 直接降低边的流量上限
    for a, b in edge_net.keys():
        vol1, cap1 = edge_net[(a, b)]
        vol2, cap2 = edge_net[(b, a)]
        min_cap = get_min_cap((a, b), len(vs))
        if vol1 < min_cap + 1000 and vol2 < min_cap + 1000:
            edge_net[(a, b)][1] = edge_net[(b, a)][1] = min_cap

    # 如果流量可以进入地面并且不造成拥堵，将多余的流量放入地面
    # 如果不能全部放入地面，将不能放入地面的部分试着分流
    # 如果不能分流，增加边的流量上限
    for a, b in edge_net.keys():
        vol, cap = edge_net[(a, b)]
        if vol > cap:
            over = vol - cap
            up_bound = 0
            if a < len(vs) and b < len(vs):
                up_bound = min([self_outs_up[a], self_ins_up[b]])
            # 放入地面
            if 0 < over < up_bound:
                edge_net[(a, b)][0] -= over
            # 分流
            else:
                over -= up_bound
                paths = get_paths(links, len(vs), a, b)
                if len(paths) == 1:
                    # 只有这一条边，无法增广，扩大边上限
                    edge_net[(a, b)][1] = get_cap((a, b), len(vs))
                else:
                    paths.pop(paths.index((a, b)))
                    for path in paths:
                        res_net = get_res_net(path, edge_net)
                        forward = [res_net[key][0] for key in res_net.keys()]
                        if min(forward) > over:
                            # 增广
                            edge_net[(a, b)][0] = cap
                            for k in range(len(path) - 1):
                                edge_net[(path[k], path[k + 1])][0] += over
                            for v in path[1: -1]:
                                all_ins[v] += over
                                all_outs[v] += over
                        else:
                            # 其他通路的剩余能力不足
                            edge_net[(a, b)][1] = get_cap((a, b), len(vs))


def reduce_circle_volumn(links, edge_net, all_ins, all_outs):
    # 降低一条环上的冗余流量
    for a, b in edge_net.keys():
        if a < len(vs) and b < len(vs):
            paths = get_paths(links, len(vs), a, b)
            if len(paths) > 1:
                for path in paths:
                    if len(path) > 2:
                        path += (a,)
                        backwards = [edge_net[key][0] for key in zip(path, path[1:])]
                        vol = min(backwards)
                        check_in = all([all_ins[i] - vol > self_ins[i] for i in path[:-1]])
                        check_out = all([all_outs[i] - vol > self_outs[i] for i in path[:-1]])
                        if vol > 0 and check_in and check_out:
                            for key in zip(path, path[1:]):
                                edge_net[key][0] -= vol
                                all_outs[key[0]] -= vol
                                all_ins[key[1]] -= vol

    # 删除流量为0的单向边
    for (a, b) in edge_net.keys():
        if edge_net[(a, b)][0] <= 0.0:
            edge_net.pop((a, b))
            links[a].remove(b)
    # 添加单向边，保证都是双向边
    for a, b in edge_net.keys():
        if (b, a) not in edge_net:
            edge_net.setdefault((b, a), [0, get_cap((b, a), len(vs))])
            links[b].add(a)

    reduce_edge_cap(links, edge_net, all_ins, all_outs)


def get_all_cuts(links, in_num):
    # 获取所有的割集
    def contract(s, t, mp, combine, parents, n):
        vis = [False] * n
        w = [0] * n
        e = []
        for i in range(n):
            e.append([])
        tmpj = n
        for i in range(n):
            maxw = -1
            for j in range(n):
                if (not vis[j]) and (not combine[j]) and maxw < w[j]:
                    maxw = w[j]
                    tmpj = j
            if t == tmpj:
                return s, t, w[t], e[t]
            vis[tmpj] = True
            s = t
            t = tmpj
            for j in range(n):
                if (not vis[j]) and (not combine[j]):
                    w[j] += mp[t][j]

            vs = [t]
            for v in range(n):
                if parents[v] == t:
                    vs.append(v)
            for v in vs:
                for j in links[v]:
                    if (vis[j] or vis[parents[j]]) and j != t and parents[j] != t:
                        if (v, j) not in e[t] and (j, v) not in e[t]:
                            e[t].append((v, j))

        return s, t, w[t], e[t]

    n = len(links)
    mp = np.zeros((n, n))
    for i in range(n):
        link = links[i]
        for j in link:
            mp[i][j] += 1
    combine = [False] * n
    parents = range(n)

    ans = []
    for i in range(n):
        s = t = -1
        s, t, minc, mine = contract(s, t, mp, combine, parents, n)
        if s > 0:
            combine[t] = True
            parents[t] = parents[s]
            for k in range(n):
                if parents[k] == t:
                    parents[k] = parents[s]
            ans.append(mine)
            for j in range(n):
                mp[s][j] += mp[t][j]
                mp[j][s] += mp[j][t]

    return sorted(ans, key=lambda x: len(x))


def find_connect_domain(links, s):
    # 求连通域
    ans = []

    def find_helper(a, ans):
        if a in ans:
            return
        ans.append(a)
        for x in links[a]:
            find_helper(x, ans)
        return ans

    return sorted(find_helper(s, ans))


def find_connect_domain_after_cut(links, cut_edges):
    # 求去除割集之后的两个连通域
    cut_links = []
    for a in range(len(links)):
        cut_links.append([])
        for b in links[a]:
            if not ((a, b) in cut_edges or (b, a) in cut_edges):
                cut_links[a].append(b)
    a, b = cut_edges[0]
    A = find_connect_domain(cut_links, a)
    B = [i for i in range(len(links)) if i not in A]
    return A, B


def add_edge_for_shunt(a, b, links, edge_net, all_ins, all_outs, in_num):
    # 加入一条新边，不断分流
    links[a].add(b)
    links[b].add(a)
    edge_net.setdefault((a, b), [0, get_min_cap((a, b), in_num)])
    edge_net.setdefault((b, a), [0, get_min_cap((a, b), in_num)])
    for x, y in [(a, b), (b, a)]:
        paths = get_paths(links, in_num, x, y)
        for path in paths:
            forward = [edge_net[(v1, v2)][0] for (v1, v2) in zip(path, path[1:])]
            in_bounds = [all_ins[i] - self_ins[i] for i in path[1:-1]]
            out_bounds = [all_outs[i] - self_outs[i] for i in path[1:-1]]
            vol = min(forward + in_bounds + out_bounds) * 0.5
            ab_bound = edge_net[(x, y)][1] - edge_net[(x, y)][0]
            if vol < ab_bound and vol > 0:
                edge_net[(x, y)][0] += vol
                for (v1, v2) in zip(path, path[1:]):
                    edge_net[(v1, v2)][0] -= vol
                for v in path[1:-1]:
                    all_ins[v] -= vol
                    all_outs[v] -= vol

    if edge_net[(a, b)][0] <= 0 and edge_net[(b, a)][0] <= 0:
        edge_net.pop((a, b))
        edge_net.pop((b, a))
        links[a].remove(b)
        links[b].remove(a)
        return False
    else:
        return True


def add_redundant_edge(links, edge_net, all_ins, all_outs, in_num):
    all_cuts = get_all_cuts(links, in_num)
    bound = len(all_cuts[0]) + 1
    bad_cut = []
    while True:
        cuts = get_all_cuts(links, in_num)
        min_cut = None
        for cut in cuts:
            if len(cut) < bound and cut not in bad_cut:
                min_cut = cut
        if min_cut is None:
            break
        A, B = find_connect_domain_after_cut(links, min_cut)
        add_edges = []
        for a in A:
            va = vs[a] if a < in_num else centers[a - in_num]
            for b in B:
                if not (a >= in_num and b >= in_num):
                    vb = vs[b] if b < in_num else centers[b - in_num]
                if (a, b) not in edge_net.keys():
                    d = get_distance(va, vb)
                    add_edges.append((a, b, d))
        is_bad = True
        for (a, b, d) in add_edges:
            if add_edge_for_shunt(a, b, links, edge_net, all_ins, all_outs, in_num):
                is_bad = False
                break
        if is_bad:
            bad_cut.append(min_cut)
    reduce_edge_cap(links, edge_net, all_ins, all_outs)


def optimize_net(links, edge_net, all_ins, all_outs):
    # 优化网络

    reduce_edge_cap(links, edge_net, all_ins, all_outs)

    reduce_circle_volumn(links, edge_net, all_ins, all_outs)

    add_redundant_edge(links, edge_net, all_ins, all_outs, len(vs))

    print check(all_ins, all_outs, edge_net, len(vs))

    transport_rate = calc_transport_rate(edge_net, len(vs))

    return links, edge_net, all_ins, all_outs, transport_rate


def get_plan(edge_net, years=8):
    # 获取8年的建设计划
    all_edges = []
    for (a, b) in edge_net.keys():
        if a < b:
            if b < len(vs):
                d = get_distance(vs[a], vs[b])
            else:
                d = get_distance(vs[a], centers[b - len(vs)])
            all_edges.append((a, b, d))
    all_edges.sort(key=lambda x: x[-1])
    center_edges = [(a, b) for (a, b, d) in all_edges if b >= len(vs)]
    inner_edges = [(a, b, d) for (a, b, d) in all_edges if b < len(vs)]
    mst = get_mst(vs, inner_edges)
    mst_edges = [(a, b) for (a, b, d) in mst]
    plan_edges = center_edges + mst_edges
    for (a, b, d) in all_edges:
        if (a, b) not in plan_edges:
            plan_edges.append((a, b))
    edge_cnt = len(all_edges)
    cap = edge_cnt / years
    plan = []
    for i in range(years):
        plan.append(plan_edges[cap * i: cap * (i + 1)])
    for i in range(years * cap - len(plan_edges), 0):
        plan[i].append(plan_edges[i])
    return plan


if __name__ == '__main__':
    links, edge_net, all_ins, all_outs, transport_rate = construct_net()
    links, edge_net, all_ins, all_outs, transport_rate = optimize_net(links, edge_net, all_ins, all_outs)
    plan = get_plan(edge_net)
    print len(edge_net)

    # 保留精度
    transport_rate = map(lambda x: round(x, 4), transport_rate)
    all_ins = map(lambda x: round(x, 3), all_ins)
    all_outs = map(lambda x: round(x, 3), all_outs)
    for key in edge_net.keys():
        vol, cap = edge_net[key]
        vol, cap = round(vol, 3), round(cap, 3)
        edge_net[key] = [vol, cap]

    # 写入csv
    with open('transport_rate1.csv', 'w') as f:
        writer = csv.writer(f)
        for rate in transport_rate:
            writer.writerow([rate])

    with open('in_volumns1.csv', 'w') as f:
        writer = csv.writer(f)
        for v in all_ins:
            writer.writerow([v])

    with open('out_volumns1.csv', 'w') as f:
        writer = csv.writer(f)
        for v in all_outs:
            writer.writerow([v])

    with open('edge_volumns1.csv', 'w') as f:
        writer = csv.writer(f)
        for (a, b) in sorted(edge_net.keys()):
            vol, cap = edge_net[(a, b)]
            writer.writerow([a, b, vol, cap])

    with open('plan1.csv', 'w') as f:
        writer = csv.writer(f)
        for p in plan:
            edges = []
            for a, b in p:
                vol, cap = edge_net[(a, b)]
                edges.append((a, b, vol, cap))
            writer.writerow(edges)
